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ABSTRACT 

We use TV-body simulations to measure mass functions in flat cosmological models 
with quintessence characterized by constant w with w = — 1, —2/3 and —1/2. The 
results are compared to the predictions of the formula proposed by Jenkins et al. 
at different redshifts, in terms of FOF masses as well as Abell masses appropriate 
for direct comparison to observations. The formula reproduces quite well the mass 
functions of simulated haloes in models with quintessence. We use the cluster mass 
function data at a number of redshifts from Carlberg et al. to constrain Slo, <J$ and 
w. The best fit is obtained in the limit w — ► 0, but none of the values of w in the 
considered range — 1 < w < can actually be excluded. However, the adopted value 
of w affects significantly the constraints in the Qq — 1% plane. Taking into account the 
dependence on w we find Qq = 0-32 ± 0.15 and cr s — 0.85to.i2 (68% c.L). Since less 
negative w push the confidence regions toward higher O and lower erg we conclude that 
relaxing the assumption of w = — 1 typically made in such comparisons may resolve 
the discrepancy between recent cluster mass function results (yielding rather low f^o 
and high a a) and most other estimates. The fact that high w values are preferred may 
however also point towards some unknown systematics in the data or the model with 
constant w being inadequate. 

Key words: methods: iV-body simulations - methods: analytical - cosmology: the- 
ory - cosmology: dark matter - galaxies: clusters: general - large-scale structure of 
Universe 



1 INTRODUCTION 

Recently, our knowledge on background cosmology has im- 
proved dramatically due to new supernovae and cosmic mi- 
crowave background data. Current observations favor a flat 
Universe with flo ~ 0.3 (see e.g. Harun-or-Rashid & Roos 
2001; Krauss 2003 and references therein) and the remaining 
contribution in the form of cosmological constant or some 
other form of dark energy. A class of models that satisfy 
these observational constraints has been proposed by Cald- 
well, Dave, & Steinhardt (1998) where the cosmological con- 
stant is replaced with an energy component characterized 
by the equation of state p/g — w 7^ — 1. The component 
can cluster on largest scales and therefore affect the mass 
power spectrum (Ma et al. 1999) and microwave background 
anisotropies (Doran et al. 2001; Balbi et al. 2001; Caldwell 
& Doran 2003; Caldwell et al. 2003). 

A considerable effort has gone into attempts to put con- 
straints on models with quintessence and presently the val- 
ues of —1 < w < —0.6 seem most feasible observationally 
(Wang et al. 2000; Huterer & Turner 2001; Jimenez 2003). 
Another direction of investigations is into the physical ba- 



sis for the existence of such component with the oldest at- 
tempts going back to Ratra & Peebles (1988). One of the 
promising models is based on so-called "tracker fields" that 
display an attractor-like behaviour causing the energy den- 
sity of quintessence to follow the radiation density in the 
radiation dominated era, but dominate over matter density 
after matter-radiation equality (Zlatev, Wang, & Steinhardt 
1999; Steinhardt, Wang, & Zlatev 1999). It is still debated, 
however, how w should depend on time, and whether its 
redshift dependence can be reliably determined observation- 
ally (Barger & Marfatia 2001; Maor, Brustein, & Steinhardt 
2001; Weller & Albrecht 2001; Jimenez 2003; Majumdar & 
Mohr 2003). 

From the gravitational instability point of view, the 
quintessence field and the cosmological constant play a very 
similar role: both can be treated as (unclustered) dark en- 
ergy components that differ by their equation of state pa- 
rameter, w. Technically, the equations governing the expan- 
sion of the Universe and the growth of density perturbations 
in the two models differ only by the value of w. Given the 
growing popularity of models with quintessence, in this pa- 
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per we generalize the description of the mass functions to 
include the effect of dark energy with constant w. 

The mass function of clusters of galaxies has been used 
extensively to estimate cosmological parameters, especially 
Qq and the rms density fluctuation as, usually yielding (due 
to degeneracy) a constraint on some combination of the two. 
It has been shown that this degeneracy can be significantly 
decreased by using the data at different redshifts (Carlberg 
et al. 1997a; Bahcall, Fan & Cen 1997; Bahcall & Bode 2003) 
due to different growth rates of density fluctuations in differ- 
ent models. The growth rate also varies among models with 
different w; one may therefore hope to distinguish between 
them by analyzing mass function data at high redshift. It 
is worth noting that the power spectra of models differing 
only by w are different only at very large scales, so the rms 
density fluctuations at scales of interest are almost identical 
when normalized to a$. 

The paper is organized as follows. In Section 2 we briefly 
summarize the properties of the cosmological model with 
quintessence, including the linear growth factor of density 
fluctuations. Section 3 presents the iV-body simulations and 
the mass function measured from their output. In Section 4 
we describe the transformation between mass functions in 
terms of FOF masses and in terms of Abell masses. Sec- 
tion 5 is devoted to the comparison of the theoretical mass 
functions to the data for clusters at different redshifts to ob- 
tain constraints on cosmological parameters. The Discussion 
follows in Section 6. 



2 THE COSMOLOGICAL MODEL 

Quintessence obeys the following equation of state relating 
its density qq and pressure pq 



PQ = W QQ, where — I < w < 0. 



(1) 



The case of w = — 1 corresponds to the usually defined cos- 
mological constant. 

The evolution of the scale factor a = R/Ro = 1/(1 + z 
(normalized to unity at present) in the quintessential Uni- 
verse is governed by the Friedmann equation 

da H 
where 



(2) 



/(«) 



l + no (I-l) +flb (^-l) 



-1/2 



(3) 



(4) 



and H is the present value of the Hubble parameter. 
The quantities with subscript here and below denote the 
present values. The parameter Q is the standard measure 
of the amount of matter in units of critical density and q 
measures the density of quintessence in the same units: 

£>crit 

The Einstein equation for acceleration d 2 a/dt 2 — 
—4nGa(p + g/3) shows that w < —1/3 is needed for the 
accelerated expansion to occur. 

Solving the equation for the conservation of energy 
d(f?Qa 3 )/da = — 3pQa 2 with condition (1), we get the fol- 
lowing evolution of the density of quintessence in the general 
case of w = w(a): 



1.0 



fi =0.3, q =0.7 



w=-l 
w=-2/3 
w=-l/2 




Figure 1. The linear growth rate of density fluctuations for Qq 
0.3, qo = 0.7 in three cases of w = —1, —2/3 and —1/2. 



QQ = QQ,o exp 



-3 In a + 3 



/' 

J a 



w(a)da 



(5) 



in agreement with Caldwell et al. (1998). For w = const, the 
case considered in this paper, the formula reduces to 



QQ = QQ,o a 



-3(1+™) 



(6) 



The evolution of f2 and q with scale factor (or equiva- 
lently redshift) is given by 



fi(a) = 
q(a) = 



O / 2 (a) 



got (a) 

n l+3w 



(7) 
(8) 



while the Hubble parameter itself evolves so that H(a) = 
H /[a f(a)}. 

The linear evolution of the matter density contrast 8 — 
8g/g is governed by equation <5 + 2(a/a)& — 4ttGqS = 0, 
where dots represent derivatives with respect to time. For 
fiat models and arbitrary w an analytical expression for D(a) 
was found by Silveira & Waga (1994, corrected for typos). 
With our notation and the normalization of D(a) = a for 
Q — 1 and q — 0, it becomes 

1 w-1 5 _ 3lu 1 - fi " 

1 , 1 , — 

3w 2w 6w 



D(a) 



■ a 2 



Fi 



2w 



(9) 



where 2F1 is a hypergeometric function. The solutions (9) 
for w — — 1, —2/3 and —1/2 are plotted in Figure 1 for the 
cosmological parameters Qq = 0.3 and go = 0.7. 



3 THE SIMULATED MASS FUNCTIONS 

In order to study the cluster mass functions in models with 
quintessence we have used a subset of the Gpc 3 dark matter 
iV-body simulations of Bode et al. (2001). Among the family 
of cosmological models they considered two are of interest 
for us here: ACDM (with w = -1, ft = 0.3, q = 0.7, 
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n h = 0.04, h = 0.67, erg = 0.9) and QCDM (with w = 
—2/3, as = 0.84, and all other parameters the same). Both 
models had the primordial spectral index n = 1. The w = 
— 1 and w — —2/3 simulations were carried out using the 
Tree-Particle-Mesh (TPM) code described in Bode, Ostriker 
& Xu (2000). 

For the purpose of this study, a new simulation with 
w — —1/2 (and <jg = 0.88, with all the remaining param- 
eters unchanged) was run using a newer, publicly available 
version of TPM code which includes a number of improve- 
ments over the earlier code (Bode & Ostriker 2003). The 
changes in the code which would affect the numerical results 
include improvements in the time stepping (adding individ- 
ual particle time steps within trees, and a stricter time step 
criterion) and domain decomposition (adding an improved 
treatment of tidal forces, and a new selection criterion for 
regions of full force resolution). The changes in time step 
only affect the innermost cores of collapsed objects, where 
the relaxation time is short (Bode & Ostriker 2003), leav- 
ing the mass function unchanged. The changes in domain 
decomposition affect the mass function, but only at the low- 
mass end: in the w = — 1 and w — —2/3 simulations, the 
mass function is complete only above 6 x 10 13 /i _1 Mq (100 
particles), while in the w — —1/2 run it extends down to 
1O 13 /i _1 M (16 particles). In this paper we concentrate on 
the high-mass end of the mass function, where these numer- 
ical differences have no impact. 

The identification of halos in the final particle distribu- 
tion has been performed with the standard friends-of-friends 
(FOF) halo-finding algorithm with the linking parameter 
b = 0.2 (particles are linked if their separation is smaller 
than b times mean interparticle distance). As an alternative, 
we have also used the HOP regrouping algorithm (Eisenstein 
& Hut 1998). This procedure employs several parameters; 
however, the resulting mass function is sensitive to only one 
of them, the density (in units of the background density) 
at the outer boundary of the halo, 5 utcr- As discussed by 
Eisenstein & Hut (1998), assuming <S ou tcr = 80 is equiva- 
lent to using the linking parameter of 6 = 0.2 in the FOF 
halo-finding algorithm. Indeed, we verified that the mass 
functions measured in the simulations are almost identical 
when using these two halo-finding schemes. 

Figure 2 shows in solid lines the mass functions mea- 
sured with FOF from iV-body simulations with w — — 1 (up- 
per panel), w — —2/3 (middle panel) and w = —1/2 (lower 
panel) at redshifts z = 0, 0.5, 1 and 2. The quantity plotted is 
the cumulative mass function (the comoving number density 
of objects of mass grater than M) N{> M) = n(M)dM, 
where n(M) is the number density of objects with mass be- 
tween M and M + AM. Jenkins et al. (2001) established 
that when the masses are identified with FOF(0.2) the sim- 
ulated mass function in different cosmologies and at differ- 
ent epochs can be well approximated by the following fitting 
formula 

where 

F(M) = 0.315 exp(-|lno- _1 + 0.61| 3 - 8 ). (11) 

In the above formulae Qh is the background density, a is the 
rms density fluctuation at top-hat smoothing scale R 




— 4 I 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 r 




R.O 14.5 15.0 15.5 




Figure 2. Comparison of the iV-body cumulative mass functions 
(solid lines) obtained with w = — 1 (upper panel), w = —2/3 
(middle panel) and w = —1/2 (lower panel) to the predictions 
of the Jenkins formula (dotted lines) at different redshifts z = 
0,0.5, 1 and 2. 
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D 2 (a) [ 
(2tt)3 J 



d 3 kP(k)W^ H {kR) 



(12) 



where D(a) is the linear growth factor given by equation 
(9), WTii(kR) is the top- hat filter in Fourier space and the 
mass is related to the smoothing scale by M = Ang^R 3 /2>. 

P(k) in equation (12) is the power spectrum of density 
fluctuations, which we assume here to be given in the form 
proposed by Ma et al. (1999) for flat models. For the present 
time (a = 1) the power spectrum is P(k) = Ck n T 2 (k) where 
n is the primordial power spectrum index (we will assume 
n = 1) and T(k) is the transfer function. For the case of 
cosmological constant (ACDM) we take the transfer function 
Ta in the form proposed by Sugiyama (1995) 



Tl(p) 



ln 2 (l + 2.34p) 



(2.34p) 2 
x [1 + 3.89p + 

fc/(r/iMpc- 1 ) and 



16.1p) 2 + (5.46p) 3 



r = 



(13) 

(6.7ip) 4 p i/2 

f2o/iexp[— fi b (l + 



where p 

V2h/n )] 

For models with quintessence the transfer function is 
Tq — TqaTa, where Tqa = Tq/Ta can be approximated by 
fits given in Ma et al. (1999). However, for w =fc — 1, Tqa 
differs from unity only at very large scales, i.e. very small 
wavenumbers. Thus if the spectra are normalized to as (the 
rms density fluctuation smoothed with R = 8ft -1 Mpc), 
then there is no need to introduce the correction from Ta to 
Tq in (12) because it does not affect the calculation of a. 

The predictions for N(> M) obtained from equations 
(10)- (13) for different w and redshifts are shown in Figure 2 
as dotted lines. One can see that the Jenkins et al. formula 
accurately reproduces the simulated mass function, espe- 
cially at z = 0. At higher redshifts the formula of Jenkins et 
al. slightly overpredicts the number density of haloes; a sim- 
ilar trend has been found recently by Reed et al. (2003). 
Therefore, we conclude that although originally designed 
on basis of other cosmological models, it can be considered 
valid also in the presence of quintessence. We therefore con- 
firm the result of Linder & Jenkins (2003) who also found a 
good agreement between the simulations and the predictions 
of Jenkins et al. (2001) formula in a different quintessence 
model. 



4 FROM FOF TO VIRIAL AND ABELL 
MASSES 

The masses of clusters of galaxies are usually measured as 
the so-called Abell masses, i.e. the masses inside the Abell 
radii of 1.5ft -1 Mpc. For the purpose of comparison with 
observations we need to transform the mass functions ex- 
pressed in terms of the FOF masses to Abell masses. This 
can be done assuming a density distribution inside the clus- 
ter. Since clusters are believed to be dominated by dark 
matter (e.g. Carlberg et al. 1997b; Lokas & Mamon 2003) 
their density distribution can be well approximated by the 
universal profile proposed by Navarro, Frenk & White (1997, 
NFW) for dark matter haloes 
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Figure 3. The concentration parameter at z = calculated from 
the model of Bullock ct al. (2001) for U = 0.3, q = 0.7 and 
w = —1, —2/3 and —1/2 with the remaining parameters as in the 
simulations, except that we kept erg = 0.9 for all cases. 

where the radius has been expressed in units of the virial 
radius r v , s — r/r v . The virial radius is defined as the dis- 
tance from the centre of the halo within which the mean 
density is A c times the critical density, p cr it- The value of the 
virial overdensity A c is estimated from the spherical collapse 
model and depends on the cosmological model. For flat mod- 
els with constant quintessence parameter w, its behaviour 
can be well approximated as (Weinberg & Kamionkowski 
2003) 

A c = 18tt 2 (1 + a6 b )Q.(z) (15) 

with 6 = {1/Qo - 1)(1 + zf w and a = 0.399 - 
1.309[(-w) ' 426 - 1], b = 0.941 - 0.205[(-w) ' 938 - 1]. 

The quantity c introduced in equation (14) is the con- 
centration parameter, c = r v /r s , where r s is the scale ra- 
dius (at which the slope of the profile is r -2 ). The function 
g(c) in equation (14) is g(c) = l/[ln(l + c) — c/(l + c)]. 
From cosmological Af-body simulations (NFW; Jing 2000; 
Jing & Suto 2000; Bullock et al. 2001), we know that c 
depends on the mass and redshift of formation of the ob- 
ject, as well as the initial power spectrum of density fluc- 
tuations. We will approximate this dependence using the 
toy model of Bullock et al. (2001) (their equations (9)-(13) 
with parameters F — 0.001 and K = 3.0 as advertised for 
masses M > 10 14 h~ 1 Mq). The predictions of the model 
at 2 = for cosmological models as in our simulations 
(but all normalized to as = 0.9) are shown in Figure 3 for 
masses 10 14 — 10 16 h^ 1 Mq. For higher redshifts we assume 
c(z) = c(z = 0)/(l + z), following Bullock et al. (2001). 
The dependence of the concentration on w is rather weak 
but there is a trend of larger c values for less negative w. A 
similar trend was observed in the properties of dark haloes 
obtained in the A r -body simulations by Klypin et al. (2003) , 
although for smaller masses. 

Since the NFW profile describes the halo properties in 
terms of the virial radius and the corresponding virial mass 
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M v — 47rr^A c p cr i t /3, we have to transform the FOF masses 
first to virial masses. As previously noted, the FOF parame- 
ter b — 0.2 corresponds to the local overdensity at the border 
of the halo 5 ou tor = 80. Therefore we will assume that the 
FOF masses are equal to masses enclosed by such isoden- 
sity contour, Mfop = Mgo- Using the NFW distribution we 
find that for our simulated ACDM model at z = we have 
M S o ~ M v , while for higher redshifts M 80 < M v . In the 
case of models with w — —2/3 we instead get Mgo > M v , 
and more so for w = —1/2. The differences between Mgo 
and M v for the models, redshifts and mass range consid- 
ered here are typically of the order of a few percent and do 
not exceed 20%. They depend somewhat on mass, and are 
due mainly to the dependence of A c on cosmology. Once 
the FOF masses are translated to virial masses we can ob- 
tain the Abell masses Ma using the NFW mass distribution 
following from equation (14) 



M A = M(s A ) = M v g(c) 



ln(l + cs A ) 



CS A 



(16) 



1 + CS A - 

where sa = ta/t^ is the Abell radius in units of the virial 
radius. 

Figure 4 shows the cumulative mass functions in terms 
of the Abell masses as measured in the simulations (solid 
lines) and as calculated with the Jenkins et al. formula (10)- 
(11) with the proper mass transformation (dotted lines). The 
agreement is very good for z = 0. In general, the analytic 
predictions match the iV-body results as long as the virial 
radius is near or greater than the Abell radius. However, at 
higher redshifts the predictions for lower masses are overesti- 
mated. Several effects might lead to this disagreement. Abell 
masses were measured in the simulations in the manner de- 
scribed in Bode et al. (2001). In this method, clusters cannot 
be closer together than 1 ft -1 Mpc; also, if two Abell radii 
overlap, particles in the overlap region are only included in 
one of the clusters (based on binding energy). These two 
factors may lead to fewer objects - lower mass objects being 
subsumed into larger ones. 

As discussed by Bode et al. (2001) higher resolution 
simulations produce somewhat steeper Abell mass functions 
which would agree better with the predictions. We must also 
keep in mind that the result was obtained with the assump- 
tion of NFW profile, while for smaller haloes this may not 
be the case. The haloes of mass 10 14 /i _1 Mq have only 160 
particles in our simulations and can hardly be expected to 
have a well defined density profile. We have also extrapo- 
lated the model for concentration of Bullock et al. (2001) 
both for quintessence and very large masses, two regimes 
where it has never been tested, while even in well studied 
models concentration shows substantial scatter. Given the 
number of approximations involved in the result we con- 
clude that the fits are satisfactory, especially in the mass 
range between the two vertical dashed lines in each panel 
of Figure 4, covered by the data, which will be used in the 
next Section. 
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5 COMPARISON WITH OBSERVATIONS 

For comparison of the theoretical mass functions with ob- 
servations we used the data for cluster mass function from 
Carlberg et al. (1997a). The data consist of seven data points 



Figure 4. Cumulative mass functions in terms of Abell masses 
measured in iV-body simulations (solid lines) and obtained from 
the predictions of the Jenkins formula with mass transformation 
(dotted lines) for w = —1, w = —2/3 and w = —1/2 and at 
different redshifts. The two vertical dashed lines in each panel 
indicate the mass range of the data used in Section 5. 
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Figure 5. The data of Carlbcrg ct al. (1997a) for H = 0.3 
versus Abell mass functions predicted for different models with 
quintessence. The filled circles mark the four data points at low 
redshift < z < 0.1. The open symbols correspond to higher red- 
shift data: 0.18 < z < 0.35 (open triangle), 0.35 < z < 0.55 (open 
square) and 0.55 < z < 0.85 (open circle). The lines show the 
models i.e. our mass functions in terms of Abell mass discussed 
in the previous Section. All models have Qo = 0.3, go = 0.7 and 
as = 0.9 but different w. w = — 1 (solid lines), w = —2/3 (dashed 
lines) and w = —1/2 (dotted lines). The three upper lines are for 
z = 0, while the three lower ones for z = 0.5. 



of N(> M) for different Abell masses M > 4 x 10 14 /i _1 M Q , 
for redshift bins in the range < z < 0.85. There are four 
data points coming from different surveys in the low redshift 
bin < 2 < 0.1 and three data points in higher redshift bins 
0.18 < z < 0.35, 0.35 < z < 0.55 and 0.55 < z < 0.85, re- 
spectively (see Table 1 of Carlberg et al. 1997a). The data 
for fio = 0.3 are shown in Figure 5 with the filled circles cor- 
responding to the low-redshift samples and open symbols to 
the higher- redshift ones. 

Together with the data we show in Figure 5 our model 
Abell mass functions discussed in the previous Section. All 
models have fio — 0.3, go = 0.7, as = 0.9 and differ only in 
the value of w. The three upper lines shown in the Figure 
are for 2 = 0, while the three lower ones for z = 0.5. We note 
that the small difference between the models with different 
w at 2 = is due only to our mapping procedure between 
the FOF and Abell masses which involves «;-dependent pa- 
rameters A c and c. (The predictions of the Jenkins et al. 
2001 formula for FOF masses applied without any correc- 
tions would be identical for the given set of cosmological 
parameters since the differences in the power spectrum for 
different to are negligible.) The differences between predicted 
Abell mass functions for different w start to be more pro- 
nounced at higher redshifts (three lower curves in Figure 5) 
due to different growth factors of density perturbations in 
these models (see Figure 1). 

We performed a standard \ 2 fit to the data points in 
terms of log N(> M) with three free parameters: fio, as and 
w. Since the data are in the form of a log N(> M) value per 



redshift bin, for the predictions we take the mean redshift of 
the bin (it is not necessary to average over redshift because 
log N(> M) is approximately linear in redshift). When con- 
sidering different fio in models we adjust the data point by 
linearly interpolating between the data given in Carlberg 
ct al. (1997a) for fio = 0.2 and fio = 1- The analysis has 
been done only for fiat models, i.e. we kept go = 1 — fio ■ We 
also assumed Hubble constant h = 0.7 and the primordial 
spectral index n = 1. 

Figure 6 shows the la, 2a and 3a probability contours 
in the a s — fio (left column), a s — w (middle column) and 
fio — w (right column) parameter planes respectively. For 
each plane the cuts through the confidence region are done 
for three values of the third parameter; the value is indi- 
cated at the corner of each panel. The contours correspond 
to Ax 2 = X 2 ~ Xmin = 3.53, 8.02, 14.2, where the minimum 
value Xmin = 10-1 is obtained for fio = 0.32, as = 0.85 and 
w -> 0. 

The contours in the as — fio plane presented in the 
left column of Figure 6 are the cuts through the confidence 
region at (from top to bottom) w = — 1, w = —0.5 and 
w — —0.2. They show a typical shape obtained in this kind 
of analyses. However, the three cuts through the confidence 
space shown in this column of Figure 6 actually move sig- 
nificantly when the assumed value of w is changed (note 
that the axes scales in the left column of the Figure are the 
same in each panel). Taking into account the dependence 
on w and the variability of the contours in the whole range 
-1 < w < 0, we find fi = 0.32 ± 0.15 and a 8 = 0.85±g;fl 
at 68% confidence level. While for w — —0.2 the best-fitting 
values of the remaining parameters are fio = 0.3, as = 0.9, 
for lower w the contours move towards lower fio and higher 
as- For w = —1 they are centered on fio ~ 0.2, as ~ 1.1. 

The middle and right columns of Figure 6 show the con- 
straints on w in two planes, as — w and fio — w respectively. 
The middle column has the cuts for fio = 0.2, fio = 0.3 and 
fio = 0.4 from the top to the bottom panel, while in the right 
column the values of the third parameter are as = 0.8, 1 and 
1.2. As can be seen in the plots, there is a strong degeneracy 
between w and any of the two remaining parameters. The 
particular shape of the contours can be understood by re- 
ferring back to Figure 1, showing the growth rate of density 
fluctuations for different w. By changing the normalization 
of the curves in Figure 1 to give the same value of D(a) at 
present (a = 1), it is easily seen that the magnitude of den- 
sity fluctuations drops faster with redshift for more negative 
w. Thus in order to reproduce a given redshift dependence 
of the cluster mass function data, for a more negative w and 
a given fio (as) a higher value of as (fio) is needed (as both 
these parameters enhance the growth of structure). 



6 DISCUSSION 

The confidence regions we obtained in Figure 6 in the as — fio 
plane are qualitatively similar to the results of the analysis 
of the same data by Carlberg et al. (1997a), which differed 
from our approach in that they used the Press & Schechter 
(1974) approximation, a power-law distribution of mass in 
clusters, and did not consider the dependence on w. While 
for w = —0.2 the best-fitting values of the remaining param- 
eters are fio = 0.3, as = 0.9 in very good agreement with 
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Figure 6. The la, 2a and 3a probability contours in the erg — Qq (left column), ag — w (middle column) and Qq —w plane (right column) 
resulting from the fits to the cluster mass function data of Carlberg et al. (1997a). The assumed value of the third parameter is given in 
the corner of each panel. 



other estimates (e.g. from CMB, see Spergel et al. 2003 or 
type la supernovae, see Tonry et al. 2003), for lower w the 
contours move towards lower f2o and higher as . At w = — 1 
they are centered on fio ~ 0.2, as ~ 1.1. 

A very similar methodology is provided by using 
Sunyaev-Zeldovich cluster surveys to probe the evolution of 
the surface density of clusters as a function of redshift (Bat- 
tye & Weller 2003 and references therein). The constraints 
derived by the Sunyaev-Zeldovich surveys are based on the 
same physical processes and models as those employed here, 



namely the mass spectrum of haloes in the quintessential 
cosmology. Battye & Weller (2003) studied the constraining 
power of Sunyaev-Zeldovich surveys by analyzing mock cat- 
alogs of future surveys. Comparing the results obtained here 
and by Battye & Weller we find that apart from the future 
PLANCK survey the Carlberg et al. (1997a) data provide 
tighter or comparable constraints to the Sunyaev-Zeldovich 
surveys. The PLANCK cluster survey will reduce the uncer- 
tainties in as and fio by roughly a factor of 2 in comparison 
with the present results. 
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Assuming the ACDM cosmological model, Bahcall et 
al. (2003) used data from the SDSS collaboration to esti- 
mate that for fio « 0.2, as ~ 0.9, a lower value than found 
here. However, the current result includes higher redshift 
data, which leads to higher estimates for as, independent of 
fio (Bahcall & Bode 2003). It seems that assuming w = — 1 
when considering cluster abundances leads to rather low val- 
ues of fio (or alternatively, high values of as); as pointed out 
by Oguri et al. (2003), these are in mild conflict with most 
other estimates. Oguri et al. (2003) suggested that decaying 
cold dark matter may resolve this discrepancy. Our results 
offer another possibility: they show that the constraints on 
fio and as from cluster mass functions are in better agree- 
ment with other estimates if the assumption of w — —1 is 
relaxed and a less negative value of w is adopted. 

We have demonstrated the existence of a strong degen- 
eracy between w and any of the two remaining parameters, 
fio and as , which is not broken in spite of using high redshift 
data: reproducing the evolution of the cluster mass function 
data requires a higher value of fio or as for more negative w. 
A similar behaviour of the confidence regions, including the 
degeneracies, in the as — w and fio — w planes was recently 
observed also by Schuecker et al. (2003, see the lower panels 
of their Fig. 3) who used the REFLEX X-ray cluster sam- 
ple. As discussed by Douspis et al. (2003) and Crooks et al. 
(2003) the degeneracies between w and other cosmological 
parameters also plague the constraints obtained from other, 
e.g. CMB, data sets. 

The best fit to the cluster mass function data is obtained 
for a surprisingly high value of w « 0, but the dependence 
on w is rather weak and no value in the entire considered 
range — 1 < w < can actually be excluded: for every value 
in this range a reasonable combination of fio and as can 
be found which places the point in la confidence region. 
It is interesting to note, however, that depending on the 
values of fio and as, we get upper or lower limit on w: for 
high fio and low as the contours tend to provide a lower 
limit on to, while for low fio and high as we have an upper 
limit. Combined with additional constraints on fio or as 
from other data sets, cluster mass functions can therefore 
prove useful in estimating the value of w. Such an analysis 
has been recently performed by Schuecker et al. (2003) who 
combined the data from X-ray clusters with the data for 
SNIa. Since the supernova data show a strong preference for 
negative w the resulting best-fitting value of w is very close 
to -1. 

With the estimates of the cosmological parameters 
presently available the analysis presented here tends to pro- 
vide a lower limit on w. For example, for the best estimates 
from WMAP (Spergel et al. 2003), fi = 0.27 and a s = 0.9, 
the constraint from our analysis is —0.5 < w < at 95% 
confidence level. The values of fio and as, however, are not 
yet known exactly so a proper combination with other data 
sets would have to be performed. Typically, other data sets 
give an upper limit on w, e.g. WMAP in combination with 
other astronomical data gives w < —0.78 (Spergel et al. 
2003) while from the SNIa data Tonry et al. (2003) obtain 
w < —0.73. This means that the combination of these limits 
with the cluster mass function data will probably result in a 
preferred range of w instead of only an upper or lower limit. 
The high values of w preferred by the cluster mass function 
data may also turn out to be in conflict with other estimates. 



This would point towards some unknown systematics in the 
data or inconsistencies in the models. Then the cosmolog- 
ical model with constant w would have to be rejected and 
replaced with one involving some time-dependence of w. 
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